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^ , Abstract 
■ 

, We present a systematic way to analyze and model systems having many 

• characteristic time-scales. The method we propose is employed for a test- 

' case of a meandering jet model manifesting chaotic tracer dispersion with 

, long time-correlations. We first choose a suitable state space partition and 

' analyze the symbolic dynamics associated to the fluid particle position. 

' In a second step we construct a stochastic process in terms of a multi-time 

. Markovian model. This corresponds to a hierarchy of random travelers 

• 1-H ' on a graph where each traveler moves at his own time scale. The results 
pH ^ are compared on the basis of statistical measures such as entropies and 
• • ■ correlation functions. 

> ■ 

%j ■ 1 Introduction 

. d ; 

Extended systems, as they are found in nature, often show complicated dy- 
namics on several characteristic scales and times. One, often cited, example is 
turbulent motion, where the dynamics stretches from an integral scale to the 
dissipation scale. In many situations, it is possible to find a low-order approxi- 
mation of the dynamics and one remains with the task to describe the temporal 
behavior of the system under consideration. A standard example is the area 
of pattern formation [Q where one ends up with systems of ordinary differen- 
tial equations, describing the evolution of a few spatial modes. It is, however, 
not always possible to perform a mode decomposition, but nevertheless spatial 
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structure is clearly present. Then it makes sense to identify states and investi- 
gate the dynamics and interaction of these states. There exists a huge amount 
of works concerning deterministic and stochastic modeling Examples 
for the above mentioned systems are found among others in climatology Q, 
oceanography ||], or ecosystems In these examples, one typically follows 
tracer particles that are put into the respective flow and performs some parti- 
tioning afterwards. We do not want to bother with this first, very subtle step 
of partitioning, rather we take the partition as given. 

On the other hand, there exist systems where one directly observes the 
temporal signal of several states, which can be coupled. A partitioning is then 
obviously not necessary and one is not bothered anymore with the problems 
arising from that point. Examples are systems of coupled oscillators, coupled 
ode's. Between the two described scenarios lie numerical simulations where 
it is often easy to identify spatial states by physical [0, ^ or mathematical 
arguments. 

Often, one considers coupled systems with one dominant time scale and 
different typical amplitudes. States with fast oscillations are usually assumed to 
have small amplitude and are neglected. This way, small oscillations that might 
well play a role in inducing finite-size perturbations and chaos in a system are 
neglected. On the other hand, it can well appear that the fast fluctuations have 
large amplitude and thus obscure the slow states dynamics. 

Chaoticity of nonlinear systems implies the loss of memory of a tracer that 
moves in the (either extended or discrete) system. It thus makes sense to use 
a statistical description. In this article, we present a way to treat systems with 
several states, which are possibly coupled and act on different time scales. We 
explain a method to filter time signals and extract information at different time 
scales along with a non-conventional, general method to construct a stochastic 
model capable to approximate the evolution of the system under observation. 
We end up with a statistical description of a multi time scale dynamics in terms 
of a Markovian cascade process. The different time scales are identified and 
the dynamics is reconstructed in the statistical sense. This is done in order to 
compare the reconstructed and original dynamics by means of complexity mea- 
sures and correlation functions to ensure that the modeling is senseful. Apart 
from the theoretical interest about how to formulate a statistical model, we can 
figure out mainly two practical uses: firstly, one can implement a statistical 
model as a building block in geophysical or other large scale simulations to get 
an idea of the fast (and possibly small) scale fluctuations and secondly, on could 
use transition times and probabilities between different states to evaluate the 
most probable traces a passive scalar (e.g., a contaminant) would have in the 
described flow (e.g., the atmosphere). 

In this article we treat an example from geophysics: the chaotic dispersion 
process generated from a meandering jet model. Apart from the interest in 
itself, we believe that this system is worth noting for the following reasons: i) it 
shows chaos and mixing on two time scales as desired by our research goal (and 
thus represents some minimal setting giving the possibility for detailed studies 
of the methods we developed), ii) it is typical for geophysical applications in) 
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we can compare results explicitly, since there exists an older article attempting 
to solve the task with conventional methods. 

The article is organized as follows: After this brief introduction, we present 
in Section |^ the model, in Section ^ the analysis method is explained, in Section 
^ the construction of the model trajectory is described. Section || contains the 
results and the comparison of constructed and original signal, finally we end 
with a short discussion and some conclusions in Section p. 



2 The model: a signal with different character- 
istic times 

Let us consider a simple but non-trivial case of multi-time dynamics in a geo- 
physical system: Lagrangian transport across a meandering jet, formerly intro- 
duced as kinematic model of the Gulf stream by Bower ||l^ and Samelson |pl| . 
The model serves our purposes, since it can be well used as the test case in 
which only two characteristic times are involved as we will see below (i.e. fast 
and slow particle transitions across the flow). In this model, the gulf stream is 
represented by a central meandering jet, flanked by gyres rotating in opposite 
directions (Fig. |^). The model is described by the stream function pH |: 

^ , , , II — Bq cos(kx) . , 

^(x,T/) = -tanh / ^ ' -.+cy, (1) 



where x and y are the spatial coordinates of a fluid particle, k is the spatial wave 
number of the meander structure, Bq controls the amplitude of the meanders, 
the term in the denominator of (|l|) defines the width of the jet and c is the 
velocity in the "far field" north and south of the (westerly) jet current. The 
fluid particle velocity components (w, v) at the point (a;, y) are: 

and v^—. (2) 

Chaos is introduced by adding an explicit time dependence to the stream 
function (|]) by replacing Bq with B(i) defined as 

B{t) = Ba + e- cos{ujt + 9) (3) 

This induces a periodic oscillation of the meander, with period T = 
amplitude e and initial phase 6, around the mean value Bq. Depending on 
the values of the perturbation parameters, uj and e, chaotic motion occurs in 
the vicinity of the separatrices. For appropriate values of the perturbation 
parameters, resonance overlap occurs, and transport between northern and 
southern gyre can take place. 

A study of Lagrangian trajectories in one elementary cell of the streamline 
pattern shows that there are two basic time scales: Tg which is the advection 
time of a particle traveling along the jet core, and Tq = 27: /lj which is the 
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Figure 1: Left: snapshot of the streamlines in the meandering jet model The 
gyres are seen as closed orbits to the north and south of the jet. The boundary 
between the jet and the gyres, as defined here, correspond to = ±0.60. Spatial 
coordinates are measured in units of the wavelength 2n/k. Right: a particle 
trajectory (periodic boundary conditions). In the corresponding gulf stream 
system the particle experiences an average drift from the left to the right. The 
tracer can be caught by a gyre (states 1 and 3) or the jet (state 2) for a longer 
time (dense tracks). 



period of the oscillation of the meander. In the following, in order to obtain 
non-dimensional equations we rescale time with Tq. The particle position is 
sampled accordingly at the end of each period. 

A coarse graining of the space can be obtained by dividing the flow do- 
main into elements, e.g., separating regions with closed from regions with open 
streamlines. The partition to be used throughout this article is defined as: 1) 
the northern gyre, 2) the jet, 3) the southern gyre (cf. Fig. Q). Assigning each 
of the partition a number yields a discrete description of the trajectory. 

An example of a particle trajectory is seen in Fig. |l|. For this example 
(and the rest of this article) the parameters of the models have been chosen as: 
meander wavelength Iq = 27r/fc = 7.5, the (mean) meander amplitude Bq = 1.2 
and c = 0.12. An example of the corresponding trajectory in terms of the 
partitions is seen in Figure ^, top line. The trajectory consists of fast oscillations 
between one of the gyres and the jet, interleaved by periods where it stays in 
one of the gyres for a longer time. The fast oscillations are basically due to 
the fact that the particle undergoes a large number of crossings back and forth 
across the boundary of two partitions. 

3 Entropic analysis and filtering 

The aim of this section is the presentation of a suitable filtering technique which 
allows to highlight and characterize the dynamical features at different charac- 
teristic times. This filtering method, developed on the basis of an exit time 
approach to entropy computation, is a way to understand how to find a good 



4 



1 



o 



0.2 



0.8 



0.6 : 



0.4 -1 












10 



20 



30 



40 



50 



60 



t/T 



Figure 2: Correlation function of the signal shown in Fig. ^, top graph. The 
long-time correlations are clearly present but cannot be captured by a usual 
Markov chain with equidistant time steps. 



stochastic model for the process. 
3.1 Encoding of the signal 

Let us consider the time signal of an observable s{t). Typically s{t) is registered 
at a certain sampling rate t~^, i.e., we consider a discrete time evolution, s{n) 
with t — n ■ Ts- A coarse-grained description of s(n) is obtained by dividing 
the range of s{n) through a partition with, say, P elements, and assigning a 
label to each element of the partition, say cti, (T2, cp. Consequently, the 
time evolution of s{n) is mapped into a symbolic dynamics cr(n), with cr(n) CE 
{ci, (72, CTp}. Let us indicate with a a generic symbol of the sequence and 
with T the lifetime of a, i.e., the time the system takes to move from the state a 
into a new one. We will use sometimes the terms residence time or exit time as 
a synonym for "lifetime" when more suitable, in the context of markov chains, 
the term waiting time is used (mainly in the literature of queuing systems). 
A sequence of length n is designated as follows: 



For example, imagine a signal s{t) — (1, 1, 1, 2, 2, 2, 2, 3, 3, 1, 1, ...) which assumes 
positive integer values at each unit time step (ts = 1); then the corresponding 
(cr, r) variables are ai = (1, 2, 3, 1, ...) and r.^ = (3, 4, 2, 2, ...), respectively. 



5, 



n 



{xi,X2, ■ ■ ■ ,Xn) , where 



(4) 
(5) 
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3.2 Estimating Entropies 

Let us briefly remind some recent results concerning the development of the 
estimation of entropy using the exit-time approach |13, To illustrate the 
basic ideas we describe the method and afterwards discuss its application to the 
problem under consideration, some details are explained in the appendix. 

Out of the set of sequences Sn we collect statistics for each subsequence Sn 
and construct the probability distribution P(5'„). Now, one can compute the 
Shannon entropy h of the signal a by the following identity: 



h* 



(6) 



where h* is the Shannon entropy of a;„, h* = lim„^oo ^r*+i ^ (^^e informa- 
tion production with n), with 

H: = -J2PiSn)\nPiSn), (7) 
and (r), the average life time, is computed as 

(r) = lim — V r, (8) 

i=l 

Denoting the entropies of the a and r sequences with h* (cr) and h* (r) , respec- 
tively, one can determine upper and lower bounds on the entropy: 

max [h{a),h{T)] < h{a, t) < h{a) + h{T) (9) 

The details of the method are described in [|l3|, Q . 



3.3 The filtering procedure 

Let us now introduce a procedure to select the temporal behavior "slower than 
a given frequency" . The filter should therefore be able to discard the many fast 
fluctuations from one state to another which occur very frequently (see figure ^ 
top line). We call this filter "killing window"; it operates in the following way. 
If the signal fluctuates from one state into another state before returning to 
the first state, and if this fluctuation lasts for a shorter time than Tp (the filter 
length of the killing window), this fluctuation is ignored. Written symbolically 
this means that: 

(cri,ri), (cr2,T2), (cti,T3) -> (cri,ri + T2 + T3) , (10) 

provided that T2 < tf- The effect of the algorithm is shown in Fig. ||, and the 
role of different lengths of killing windows is shown in Fig. ^. This filter has the 
effect of "killing" all fluctuations smaller than Tp. Note that is also effectively 
prolongs the slower oscillations, c.f. Fig. 0. 
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Figure 3: An illustration of how the killing window algorithm changes the orig- 
inal signal (top) to a filtered signal (bottom). In this example the signal is only 
alternating between two states, 1 and 2. 




Figure 4: Comparison of signals filtered with different killing windows, tf = 
0, 2, 4, 8, 12, 20, from top to bottom. It is evident how the fast fluctuations are 
removed by the killing window. For a very long killing window, tf = 20, also 
the residence in the jet is killed, and the model has effectively only two states. 
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Figure 5: Calculation of the entropies, a) the entropy for Tp — 2 together with 
the upper and lower bounds form (^) . b) Comparison of the entropy for varying 
widths of the killing window. The entropies are normalized by adding a constant 
such that hi = 1. 



With this filtering procedure applied to the original signal a;„, we obtain a 
filtered sequence x'!^'' , for which we can compute the entropy h{Tp) using (^. 
We used for simplicity a filter that acts backwards in time. It is however not 
difficult to redefine the filter such as to be invariant under time reversal, the 
result remains the same within an error of the fast time scale. The observed 
signal varies over two very distant characteristic time scales, Tjast and TsIow, 
respectively. Thus, information at the slower time scale can be extracted if 
Tsiow ^ Tp Tfast, bccausc duc to the second inequality almost all fast oscilla- 
tions are killed. 

In Fig. ^ the convergence of the entropy, at varying the block size, for a 
killing window of length 2 is shown, together with the upper and lower bounds 
from (|). First of aU we see that the entropy lies nicely in between the upper 
and the lower bounds. For a killing window of length 2 this is not surprising, 
as most of the transitions are fast oscillations between jet and gyre. However, 
this entropy converges after n — 2, which means that one needs a second order 
Markovian process in order to properly describe this process. Let us stress that 
now the state of the system is given by a and r, therefore a Markov model of 
order m corresponds, in the original space of the symbols cr, to a Markov model 
of order m ■ (t/ts). The entropy based only on the residence times, /i(r), is 
larger, and thus constitutes the effective lower bound for the total entropy. In 
Fig. ||b this is shown for the bare signal, and for killing windows of a length from 
2 to 20. This entropy based on the residence times has converged immediately, 
which means that the residence times are basically uncorrelated. This is a sign 
for an exponential distribution to be confirmed below. Now, we have to decide 
the filter length for which the signal is described by a process on one time scale 
only. The killing window clearly destroys information. If the loss of information 
is proportional to the filter length, obviously the process is homogeneous (with 
respect to information loss) within the time scale we are considering. This 
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Figure 6: Information loss with filter length. From Tp — S we have a con- 
stant rate of information loss and thus the difference between the different filter 
lengths is zero. This in turn justifies the description of the process by one scale 
only. 



is analogous to the information production rate which has to converge to a 
constant to find the entropy of a chaotic system. In Fig. ^ we display 

the entropy loss depending on the filter length. For = 8 we have certainly 
reached convergence and thus the process can be described at the level of one 
single scale (namely the slow one). For completeness, we remark that almost all 
information is destroyed if the filter length exceeds the typical slow time scale, 

4 Construction of the multiple time Markovian 
Model 

In this section we will demonstrate how the analysis of the signal can be ap- 
plied to construct a Markovian model, which can generate stochastic symbolic 
sequences with the same statistics as the original signal. The model should be 
able to reproduce both the slow and fast transitions. This is accomplished by 
first reproducing the slow scale signal, and then nesting the fast transitions into 
that one. We thus define a "horizontal" and a "vertical" dependence, the first 
describes the process on a specific time scale (slow or fast), the latter models 
the interdependencies between the different time scales. Both are assumed to be 
Markovian of first order, as our preceding analysis has shown, a generalization 
is straightforward. In the following, for the sake of clarity, we will use capital 
letters for slow scale and small ones for the fast scale. 

4.1 Slow transitions 

As determined by the Entropic analysis, the slow signal is the one filtered by a 
killing window of length > 8. The construction rules are the following: 
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i) Given an entering state a/ the process selects according to the transition 
probabihty matrix Wjj its next state aj. 

ii) Given that selection of the transition IJ, the time spent in / before jumping 
to J is given by the random variable r/, which in turn possesses the distribution 
Piri). 

Both needed informations are obtained from the filtered signal. By the given 
rules, we recognize that the construction concentrates on the transitions between 
states, rather than on the state itself. The lifetime of the states depends only 
on this transition, we see that the transition matrix does not include self-loops, 
like I ^ I and thus has zeroes on the diagonal. In the case of independent 
lifetimes, i.e., exponential distribution, one can interchange steps i) and ii) and 
first select a state. Numerical examples for the matrices and residence times are 
given Section ||. 

The reconstruction corresponds to a random traveler on a graph, which stays 
in each node with a certain node-characteristic time and then jumps over an 
edge to the next node according to a prescribed hopping probability. 

4.2 Nesting of the fast scale 

Up to now, the construction is rather standard. A new thing we introduce is the 
vertical connection between scales. The dynamics of the fast scale shall depend 
on the slow one, but apart from that be a first order Markovian process as 
described above by the rules i) and ii). The two quantities we have to consider 
are the pdf of the lifetimes p{Ti) and the transition probabilities Wij. Both can 
be quite different from the slow processes quantities. 

One obvious condition on the lifetimes is, that the duration of the fast states 
time intervals should not exceed the slow ones. This is fulfilled, since we an- 
ticipated a clear separation of time scales. In general, the transitions will be 
hit with an uncertainty of (r/ost) — tf- This in accordance with the filtering 
analysis, which is only accurate up to the filter length. The lifetime distribution 
function p can be taken without any restrictions from the analysis, since we 
anticipate the times to be independent from each other (as confirmed for our 
example, see next section). The transition probabilities Wij, however, have to 
be conditioned on the slow scale transitions and thus we have to write them for 
a complete description as Wij{IJ). Results are given in the next section. 

We want to remark, again, that a first order in the space a, r corresponds 
to a very high order in the space a with conventional time sampling. For our 
illustration by the traveler on a graph we imagine the traveler moving with his 
own, slow speed, being accompanied by his dog which moves forth and back 
between the present node of the traveler and some other nodes he's allowed to 
visit. The dog rests only for short times in the nodes, and thus effectively moves 
at a fast time scale. 

In conclusion, we have to find for each of the transitions IJ the according 
fast transition probability Wij{IJ). The fast process evolves then with the rules 
from Wij {IJ) and the lifetime statistics p{t) within the window between two 
transitions. 
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Figure 7: The pdf of the residence times in the jet and the gyres, calculated 
with a kiUing window of length 2 (left) and 12 (right). 



In a generalized model, N different time scales are involved, then we order 
the processes from slowest to fastest time and construct the hierarchy by nesting 
the scale n into scale n — 1. This generalization is straightforward and will be 
reported elsewhere. Another generalization concerns Markov processes of higher 
order. This requires accurate bookkeeping of the dependencies but in principle 
is straightforward, too. 



5 Results 

In this section we will present how the technique described above compares to 
the real signal. The slow scale transition matrix, we obtain for our model 



Wij = 0.5 0.5 . (11) 




Note that due the symmetries of the meandering jet problem, this matrix could 
in principle be reduced to only its upper half. In general, this is however not 
the case. Since we have a very simple system with strong restrictions the fast 
transition matrices look simple, e.g., for 10^ (12): 

/O 1 0\ 

w,j{\2)^ 10 0. (12) 
\ / 

The residence times for the slow and fast scale are displayed in Fig. |^. Clearly 
they are very close to be exponential. This is a sign of a non-correlated process, 
which has also been found by the entropic analysis. From a fit to an exponential 
distribution, we determine the mean residence time for the slow time scale to 
be {T)siow,gyre = ^4.5 and {T)gig^ .j^t = 90.9, for the fast time scale we obtain 
{T)fasuavre = 5. 21 and (T)^„,t^^gt = 5.83, respectively. 

For our example of a coarse grained spatial system, the trajectory, basically 
slowly drifts forth and back between north and south gyre, crossing the jet. 
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During this evolution, fluctuations occur about the boundaries of the states. 
Take, e.g., a trajectory starting at the northern gyre. Having made the transition 
to the jet on the slow scale, still, fast fluctuations will occur between the gyre 
and the jet, thus masking the signal as in Fig. ^ If the trajectory evolves 
towards the southern gyre, the fluctuations will happen only between southern 
gyre and jet. 

This does not at all conflict with a probabilistic interpretation: The signal 
chooses randomly, according to the transition matrix and lifetime distribution 
a new state and lifetime. Information about history and details of the space 
evolution are forgotten (up to 1st order). This is the natural result of the coarse 
graining. The fast fluctuations are again random and follow their description 
with matrices and lifetimes. In our case, we find the trivial transition matrices, 
like in Eq. ( |l2[ ) due to the simplicity of the meandering jet system. In general, 
however, more complicated scenarios will appear. 

In Fig. ||b we show the signal for Tf — 2 and a model trajectory for the same 
parameters, as constructed by our method. In contrast to an older work |^ the 
signals resemble each other a lot, for comparison these results are reproduced 
in Fig. 

Now we present some more quantitative comparisons. The first quantity to 
check is the information contained in the signals, it should be approximately 
the same, especially if we consider that this has been the main tool to find the 
models order. In Fig. ^ the comparison for tj = 2 and tj = 12 is shown - we 
find an excellent agreement. 

Finally we check for the correlation function as the standard tool to inves- 
tigate dependencies within the signal. The result is displayed in Fig. |l^. The 
general behavior for the correlation functions is the same. It seems that the 
reconstructed ones slightly underestimate the correlations in the signal. There 
are two mechanisms leading to deviations of the construction from the original: 
firstly, the result depends sensitively on the correct estimation of the time scales 
involved in the original signal and secondly, the approximation of the pdf for 
the residence times by an exponential is nothing exact, Nevertheless, the coin- 
cidence is remarkably good, taking into account the possible error sources. For 
real-world applications, one should naturally calculate the pdf numerically and 
use the stored values. 

6 Discussion and conclusion 

In this paper we present a method to analyze and reconstruct multi-time-scale 
signals, which originate either from extended (after partitioning) or discrete 
systems. In our example, we use a simple partition of an extended system. 
To extract different time scales of the process, we apply the killing window 
filtering technique to the signal. With the filtered signals we use exit-time 
entropies to determine the Markovian order for the stochastic description of 
the process under consideration. With the information from the previous steps, 
we construct the multi time scale signal. The reconstructed signal shows good 
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Figure 8: Comparison of the reconstructed signal with the true one. a) conven- 
tional construction, using fixed sampling time (cf. [||). Upper part: original, 
lower part: construction. The slow scale transitions are not recovered and the 
signals differ a lot. b) signal as reconstructed by the Markovian cascade. Upper 
part: original, lower part: reconstruction. It is not obvious at once, which signal 
is the original one. 
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Figure 10: Left: correlation function for the filtered signal with Tp = 12 and 
the corresponding reconstructed signal. Right: the same for tf = 2. The rapid 
decay of the correlations is estimated very well by the model data. 
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statistical coincidence with the original data. 

We did not investigate in detail the question about the amount of data 
needed for our method. This might vary from case to case, according to the 
conditions that are needed for the respective application. But error estimates 
are standard in this case and we do not expect novelties or surprises. The 
procedure is consistent as a whole, there are no ad-hoc assumptions necessary 
and everything is based purely on observations. 

More general scenarios with more than two time scales involved are straight- 
forward generalizations of the presented scheme, dependencies at a deeper level 
can be included easily, e.g., dependence of a fast time scale on the two preceding 



slower ones. This concept resembles the slaving idea |16, |17| which proved to 
be a powerful tool in pattern formation. 

Applications are, to our opinion, mainly in the modeling and forecasting of 
signals. One can well imagine such a stochastic model as a minimal building 
block within a GCM or other simulation tools, which cannot resolve fine-scale 
temporal motion due to numerical restrictions. In this way, short-time fluctua- 
tions can be modeled by data analysis with signals from experimental data. 

As a stand-alone tool, one can use the obtained information to predict the 
most probable track a tracer takes within a certain time, in this case the recon- 
struction method is only good for supervising the quality of the reconstruction. 

It is hoped that these tools can be applied to to real data. We think es- 
pecially at data of the trajectories of Lagrangian tracers from oceanographic 
experiments. There, we have to occupy with finding the right partitioning of 
space. One can try to use ideas of dynamical systems theory, as e.g. described 
in 1^. These concepts however stand and fall with the amount and quality 
of data. If more information about the system is available to yield additional 
reasoning for specific partitions this should be used as input for a meaningful 
coarse-graining. 
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8 Appendix A 



The calculation of the block entropy in Eq. || is valid for countable symbols. 
In this case, however the symbols are not countable, as each symbol involves a 
continuous variable; the residence time. The block entropy should therefore be 
calculated as 

^" = E / log(p('^") ' (13) 

where the sum extends over all possible permutations of n'th order. For practical 
purposes, the residence times have to be binned, the usual choice being a binning 
with resolution r^. Then the block entropy turns, after renumbering all the 
possible words, into the sum H"{Tr) = J^Pi^^) logp(*S'"), with summation 
over all possible words now in in the now doubly discrete space. The entropy is 
found as 

/i(t^)" = lim iJ(T^)"+^ - HiTrY' . (14) 

n — >oo 

The limit of infinite time-resolution gives us the entropy per exit, i.e.: 

/i" = lim h"{Tr) . (15) 

In Ref. [^3[ Q more details are given, including the topics of continuous 
signals and several rigorous bounds for the entropy together with some applica- 
tions. 
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